
## ---------------- ##
## Helper functions ##
## ---------------- ##

## function to pool standard errors based on Rublin's rule
rubin.pool.se <- function(se, est) {
  var_w <- mean(se^2)
  var_b <- sum((est - mean(est))^2)/7
  var_total <- var_w + var_b + var_b/7
  pool.se <- sqrt(var_total)
  pool.se
}

## Processing jackknifing models
make.texreg.jackknifing.models <- function(fitted.model.list = fitted.models.Table1, 
                                           jackknifing.ergm.fit = jackknifing.ergm,
                                           jackknifing.ergm.gof.object = jackknifing.ergm.gof) 
  { ## make texreg object from original models
  jackknifing.models <- lapply(sapply(fitted.model.list, "[", 1), texreg::extract.ergm)
  
  ## loop over each issue area in selection
  for (x in 1:length(selection)) {
    ## make temp object for each issue area
    jackknifing.issue.temp <- 
      jackknifing.ergm.fit[topic == selection[x], .(est = mean(est), 
                                         se = rubin.pool.se(se = se, est = est)), 
                           by = vars] %>% .[, .(vars, est, se, p.val = 2*pnorm(-abs(est/se)))] %>%
      .[order(match(vars, names(coef(fitted.model.list[[1]]$fit)))), ]
    
    ## replace the texreg object attributes from temp object
    jackknifing.models[[x]]@coef <- jackknifing.issue.temp[, est]
    jackknifing.models[[x]]@se <- jackknifing.issue.temp[, se]
    jackknifing.models[[x]]@pvalues <- jackknifing.issue.temp[, p.val]
    jackknifing.models[[x]]@gof <- as.numeric(
      jackknifing.ergm.gof.object[topic == selection[x], .(mean(aic), mean(bic), mean(lik))])
  }
  
  names(jackknifing.models) <- selection
  ## return the replaced texreg object
  jackknifing.models
}

## Custom GOF function for valued ERGM
gof.valued.ergm <- function(fit) {
    estsim1 <- simulate(fit,
                            monitor = ~ cyclicalweights("min", "max", "min") + 
                              transitiveweights("min", "max", "min") + 
                              mutual(form = "min", threshold = 0) + 
                              nonzero, nsim = 1000,  output = "stats")
    ## normalize the statistics based on observed mean from simulations
    estsim <- apply(estsim1[, 9:12], 2, scale) 
    estsim <- reshape::melt(estsim) ## convert to long data format for ggplot2
    estsim <- estsim[, -1]
    colnames(estsim) <- c("terms", "y")
    #This will shows you descriptive statistics in the actual observed network
    estobs <- summary(
      g ~ cyclicalweights("min", "max", "min") + 
        transitiveweights("min", "max", "min") + 
        mutual(form = "min", threshold = 0) + 
        nonzero, response = "edge_count")
    estobs <- (estobs - apply(estsim1[, 9:12], 2, mean))/apply(estsim1[, 9:12], 2, sd)
    estobs <- data.frame(terms = names(estobs), y = estobs)
    rownames(estobs) <- NULL

    out <- list(estsim = estsim, estobs = estobs)
    out
}
    
    
## --------------- ##    
## DATA PROCESSING ##
## --------------- ##

dat <- haven::read_stata("ZA5864_v1-0-0_Stata14.dta")
dat <- as.data.frame(dat)

## Include organizations
## SPÖ = 11, ÖVP = 12, FPÖ = 13, GRÜNE = 14, BZÖ = 15, STRONACH = 16 
orgas <- 11:16
dat <- dat[dat$v25b_l1 %in% orgas,]
dat <- dat[dat$v30b_l1 %in% orgas,]

## Only newspapers
dat <- dat[dat$v3_l1 == 3000,] # Tagespresse

## extract topics as two digit codes 
dat$v27_l1 <- as.numeric(str_extract(dat$v27_l1, "\\d{2}"))

## create "media" data frame
media <- data.frame(
  subject = dat$v25b_l1,
  object = dat$v30b_l1,
  valence = dat$v31,
  topic = dat$v27_l1,
  media = dat$v3_l3,
  stringsAsFactors = F
)

unique_actors <- unique(c(media$subject, media$object))
unique_topics <- unique(media$topic)


#-----------------------
# Press releases
#-----------------------

## Pure organizations are a 1 in the data set

# subject actors: v6, v8
# subject actor organizations: v7, v9
# object actors: v18, v23, v28
# object actors organizations: v19, v24, v29

ots_dat <- haven::read_stata("ZA5861_v1-0-0.dta")
ots_dat <- as.data.frame(ots_dat)

## Convert to integer, extract organization code
actors <- c(7, 9, 19, 24, 29)

for(i in actors){
  ots_dat[, str_c("v", i)] <- as.integer(ots_dat[, str_c("v", i)])
  ots_dat[, str_c("v", i)] <- str_extract(ots_dat[, str_c("v", i)], "\\d{2}")
}

## Topic
ots_dat$v11 <- str_extract(ots_dat$v11, "\\d{2}")

## Object actors 20, 25, 30 valence

## Set up edge list
ots <- data.frame(
  subject = c(
    rep(ots_dat$v7, 3),
    rep(ots_dat$v9, 3)
  ),
  object = rep(c(
    ots_dat$v19,
    ots_dat$v24,
    ots_dat$v29
  ), 2),
  valence = rep(c(
    ots_dat$v20,
    ots_dat$v25,
    ots_dat$v30
  ), 2),
  topic = rep(ots_dat$v11, 6),
  stringsAsFactors = F
)

ots$subject <- as.numeric(ots$subject)
ots$object <- as.numeric(ots$object)
ots$topic <- as.numeric(ots$topic)

ots <- ots[complete.cases(ots),]

## Discard cases that are not part of media_tmp
ots <- ots[ots$subject %in% unique_actors,]
ots <- ots[ots$object %in% unique_actors,]
ots <- ots[ots$topic %in% unique_topics,]


#-----------------------
# Shorten and relabel topics
#-----------------------

ots_old <- ots
media_old <- media

## Negative valence
ots <- ots[ots$valence == -1,]
media <- media[media$valence == -1,]

## Convert to factor
for(i in 1:2){
  ots[,i] <- as.factor(ots[,i])
  media[,i] <- as.factor(media[,i])
}

topic_selection <- list(
  c("10", "Economy"),
  c("11", "Welfare"),
  c("12", "Budget"),
  c("13", "Education"),
  c("14", "Security"),
  c("18", "Infrastructure"),
  c("19", "Society"),
  c("20", "Environment"),
  c("21", "Institutions"),
  c("22", "Immigration")
)

unique_topics_short <- unlist(topic_selection)[c(T, F)]
topics_names <- unlist(topic_selection)[c(F, T)]

ots <- ots[ots$topic %in% unique_topics_short,]
media <- media[media$topic %in% unique_topics_short,]

for(i in 1:length(topic_selection)){
  ots$topic[ots$topic == topic_selection[[i]][1]] <- topic_selection[[i]][2]
  media$topic[media$topic == topic_selection[[i]][1]] <- topic_selection[[i]][2]
}

#-----------------------
# Manifesto-Data
#-----------------------

manifestos <- haven::read_stata("ZA6877_v1-1-0.dta")
manifestos <- as.data.frame(manifestos)

## Actors
manifestos$v05 <- str_extract(manifestos$v05, "\\d{2}")
manifestos <- manifestos[manifestos$v05 %in% unique_actors,]

## Topics
manifestos$v10 <- str_extract(manifestos$v10, "\\d{2}")
manifestos <- manifestos[manifestos$v10 %in% unique_topics_short,]

issue_emphasis <- as.matrix(prop.table(table(manifestos$v05, manifestos$v10), 1))
colnames(issue_emphasis) <- topics_names

## manifesto data coding schemes
## v05 subject actor organizations (see V19 for party name type)
## v10 issues (see V10a for detailed categories)
## v14 object actor organizations

#-----------------------
# Assign names
#-----------------------

actors <- list(
  c(11, "SPÖ"),
  c(12, "ÖVP"),
  c(13, "FPÖ"),
  c(14, "GRÜNE"),
  c(15, "BZÖ"),
  c(16, "STRONACH")
)

levels(media$subject) <- 
  levels(media$object) <- 
  levels(ots$subject) <- 
  levels(ots$object) <- unlist(actors)[c(F, T)]

rownames(issue_emphasis) <- unlist(actors)[c(F, T)]


#-----------------------
# Pre-Post
#-----------------------

## Impossible to estimate field specific position
## for two issues only with one party missing and one dimension unclear
prepost <- haven::read_stata("ZA5859_de_v2-0-1.dta")

## mean left-right placements
lr <- prepost[, str_detect(colnames(prepost), "w1_q11x")]
lr[lr > 20] <- NA
lr <- apply(lr, 2, mean, na.rm = T)
names(lr) <- c("SPÖ", "ÖVP", "FPÖ", "BZÖ", "GRÜNE", "STRONACH")
lr <- lr[c("SPÖ", "ÖVP", "FPÖ", "GRÜNE", "BZÖ", "STRONACH")]

## most important question facing this country
mip <- c(prepost$w1_q5x1, prepost$w1_q6x1)
mip <- str_extract(mip, "\\d{2}")
mip[mip > 30] <- NA

## competence issue ownership
competence <- c(prepost$w1_q7, prepost$w1_q8)
competence[competence == 7 | competence == 8 | competence == 10 | competence == 11 | competence == 12 | competence == 88 | competence == 99] <- NA
competence[competence == 4] <- 3
competence[competence == 1] <- "SPÖ"
competence[competence == 2] <- "ÖVP"
competence[competence == 3] <- "FPÖ"
competence[competence == 5] <- "BZÖ"
competence[competence == 6] <- "GRÜNE"
competence[competence == 9] <- "STRONACH"

competence <- table(competence, mip)
colnames(competence) <- c(
  "Economy",
  "Welfare",
  "Budget",
  "Education",
  "Security",
  "Army",
  "Foreign",
  "Europe",
  "Infrastructure",
  "Society",
  "Environment",
  "Institutions",
  "Immigration",
  "Government"
)

competence <- competence[rownames(issue_emphasis),]
competence <- competence[, colnames(competence) %in% colnames(issue_emphasis)]

issue_emphasis <- issue_emphasis[, colnames(issue_emphasis) %in% colnames(competence)]
rownames(issue_emphasis) <- rownames(competence)
competence_prop <- prop.table(competence, margin = 2)
max <- apply(competence_prop, 2, max)

competence_bin <- competence_prop
for (i in 1:ncol(competence_bin)) {
  competence_bin[,i] <- ifelse(competence_bin[,i] == max[i], 1, 0)
}
